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Abstract. In this paper, we do a complete classification of valence-bond crystals 
(VBCs) on the kagome lattice based on general arguments of symmetry only and thus 
identify many new VBCs for different unit cell sizes. For the spin-1/2 Heisenberg 
antiferromagnet, we study the relative energetics of competing gapless spin liquids 
(SLs) and VBC phases within the class of Gutzwiller-projected fermionic wave 
functions using variational Monte Carlo techniques, hence implementing exactly the 
Q constraint of one fermion per site. By using a state-of-the-art optimization method, we 

conclusively show that the U(l) Dirac SL is remarkably stable towards dimerizing into 

fvj all 6-, 12- and 36-site unit cell VBCs. This stability is also preserved on addition 

^- of a next-nearest-neighbor super-exchange coupling of both antifcrromagnctic and 

ferromagnetic (FM) type. However, we find that a 36-site unit cell VBC is stabilized 
on addition of a very small next-nearest-neighbor FM super-exchange coupling, i.e. 
(•<-) \J 2 \ ~ 0.045, and this VBC is the same in terms of space-group symmetry as that 

*y^ obtained in an effective quantum dimcr model study. It breaks reflection symmetry, 

C^) has a nontrivial flux pattern and is a strong dimcrization of the uniform RVB SL. 

(N 



(N 



% 



PACS numbers: 71.20.-b, 75.10.Jm, 75.10.Kt, 75.30.Et, 75.40.Mg, 75.50.Ee 



Submitted to: New J. Phys. 



Valence- bond crystals in the kagome spin-1/2 Heisenberg antiferromagnet ... 2 

1. Introduction 

For many decades, physicists have been actively searching for playgrounds that are 
'hot' enough to melt magnetic freezing at temperatures well below the characteristic 
interaction energy scales in the system. This melting being fueled by quantum 
fluctuations leads to stabilization of exotic quantum paramagnetic phases of matter [1]. 
Representatives of such phases are spin liquids (SLs) and valence-bond crystals (VBCs); 
the former preserve lattice symmetries and the latter break them, according to a 
generally accepted definition. Long before any experimental hints, theoreticians such as 
Pomeranchuk already conjectured the existence of SLs [2], which were later advocated 
by Anderson to be possible appropriate ground states for the spin-1/2 Heisenberg 
antiferromagnet [3, 4]. On the experimental side, the drought in the search for SLs 
ended with the discovery of Herbertsmithite (ZnCu3(OH) 6 Cl2), a compound with perfect 
kagome lattice geometry, belonging to the paratacamite family [5, 6, 7, 8, 9, 10, 11, 
12, 13]. In it, the combination of low spin value (S = 1/2), low-dimensionality 
(d = 2) and coordination number (z = 4) and frustrating nearest-neighbor (NN) 
antiferromagnetic (AF) super-exchange interactions on a non-bipartite lattice leads to 
amplification of quantum fluctuations that stabilize a quantum paramagnet. Indeed, 
all experimental probes on Herbertsmithite point to a SL behavior down to 20 mK 
(~ J/10 4 ), which was established on the magnesium version of Herbertsmithite 
(i.e. MgCu3(OH) 6 Cl2) [14, 15, 16]. Furthermore, Raman spectroscopic studies on 
Herbertsmithite hint at a gapless (algebraic) SL [17]. 

On the theoretical side however, the nature of the ground state of the NN spin- 
1/2 quantum Heisenberg antiferromagnet (QHAF) on the kagome lattice is still elusive 
and intensely debated. Exact diagonalization studies have revealed a magnetically 
disordered ground state and a huge number of singlet excitations below the triplet 
gap [18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31]. Using approximate numerical 
techniques various claims as to the nature of the ground state have been made. These 
have included, among SL phases, a gapless (algebraic) U(l) Dirac SL using projected 
fermionic variational Monte Carlo [32, 33, 34, 35, 36], a gapped Z 2 SL [37, 38, 39] 
using density matrix renormalization group (DMRG) [40, 41] and a chiral topological 
SL using Schwinger boson mean field theory [42]. Among the VBC phases, the 
proposals have included a 36-site unit cell VBC [43, 44] numerically studied using series 
expansion [45, 46], multi-scale entanglement renormalization ansatz (MERA) [47] and 
also using quantum dimer models (QDM) [48, 49, 50]. Furthermore, VBCs with smaller 
unit cells of 6 sites [51], 12 sites [52, 53, 54] and 18 sites [43] were also argued to be viable 
ground states of the spin-1/2 QHAF. A more recent generalized QDM study found a 
new (possibly chiral) VBC of 12-site unit cell to be competing with the 36-site unit cell 
VBC. It also established an extensive quasi-degeneracy of the ground state manifold of 
the kagome S = 1/2 QHAF with a stiff competition between several phases [55]. 

In this work, we will study these non-magnetic phases within a Schwinger fermion 
formulation of the spin model. Within this approach, the projected gapless (algebraic) 
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U(l) Dirac SL has the best variational energy [32]; despite being a marginally stable 
phase, it was argued in [33] to be stable against a certain class of perturbations. Explicit 
numerical calculations using projected wave functions have in fact shown it to be stable 
(locally and globally) w.r.t. dimerizing into all known VBC perturbations [32, 34, 35]. 
Furthermore, it was shown that within this class of Gutzwiller projected wave functions, 
all the fully symmetric gapped Z2 SLs have a higher energy compared to the U(l) Dirac 
SL [36]. Similar conclusions were also reached within the Schwinger boson approach to 
the spin model [56, 57]. Note that a simple tensor network (PEPS) representation of 
such a projected bosonic RVB ansatz can be constructed and has been studied in [58]. 
In this paper, in section 2 we first perform a systematic symmetry classification of 
VBC patterns on the kagome lattice and thus identify and enumerate many new VBCs, 
independent of the formalism used to study them. In section 3, we address the question 
of relative energetics of SL and VBC phases. In particular, in section 3.1.1 we show 
that the U(l) Dirac SL is remarkably stable w.r.t. dimerizing into any of these new 
VBCs. This stability is also preserved upon addition of a finite next-nearest-neighbor 
(NNN) super-exchange coupling of both AF and ferromagnetic (FM) type. Such a NNN 
coupling might be a possible perturbation in Herbertsmithite. In section 3.1.2, we show 
that a broken symmetry phase is stabilized on addition of a small NNN FM coupling, 
which is consistent with the findings in [59]. This VBC has a 36-site unit cell with a non- 
trivial flux pattern threading its plaquettes and it is found to be a strong dimerization of 
another competing U(l) gapless SL, the so-called uniform RVB SL [32]. This 36-site unit 
cell VBC has a lower symmetry as compared to that studied in our previous work [35] 
and has precisely the same symmetry as that identified in QDM studies [49, 50, 55]. 
Thus, here we mainly establish the stability of the U(l) Dirac SL w.r.t. an extremely 
large class of potential VBC instabilities and detect a non-trivial 36-site unit cell VBC 
instability of the uniform RVB SL which is stabilized on addition of a very weak NNN 
FM super-exchange coupling to the Hamiltonian. 

1.1. The model, wave functions and the numerical technique 

The Hamiltonian for the spin-1/2 quantum Heisenberg J1 — J2 model is 

m <(«» 

where (ij) and ((ij)) denote sums over NN and NNN pairs of sites, respectively. The 
Sj are spin-1/2 operators at each site i. In the following, we will consider J% > (AF) 
and both FM and AF super-exchange J2; all energies will be given in units of J\. 

The physical variational wave functions are defined by projecting noncorrelated 
fermionic states: 

|*VMc(Xij,A#, fi,C)) ='P G \^MF(Xij,^ij,I^X)}, (2) 

where Vq, = FL(1 — ^j,t n i,4,) is the full Gutzwiller projector enforcing the one fermion 
per site constraint. Here, |\I / mf(Xjj? A»jj /x, ()) is the ground state of a mean-field 
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(a) U(l) Dirac spin liquid Ansatz 



(b) CVBC 



Figure 1: (a) The U(l) Dirac SL ansatz given up to second NN bonds. The unit cell 
has to be doubled to accommodate the 7r-flux. The black (gray) bonds denote first NN 
real hopping (second NN real hopping) terms. The solid (dashed) bonds denote positive 
(negative) hoppings. (b) The columnar VBC has no (point group) symmetries at all, 
hence all its 12 bonds are different, which are thus marked with different colors and line 
styles. Consequently, its symmetry (point) group is the identity E. 



Hamiltonian constructed out of Schwinger fermions and containing hopping, chemical 
potential and singlet pairing terms: 

H M F = ^2(Xij + ^ij)4,a C j,a + XX( Ai J + ^') C lt C L + h - C '} ' ( 3 ) 

i,3,a i,j 

where Xij — X*ji an d A^- = A^. Besides the chemical potential /i, we will also consider 
real and imaginary components of on-site pairing, which are absorbed in (. 

The SL phases are characterized by different patterns of distribution of underlying 
SU(2) gauge fluxes through the plaquettes which are implemented by a certain 
distribution of the phases of Xij an d Ay on the lattice links. Since in a SL state 
\Xij\ 2 + |Ajj| 2 is constant for each geometrical distance, a complete specification of a 
SL state up to nth NN amounts to specifying, in addition to the SU(2) fluxes, the 
optimized magnitude of hopping and pairing parameters at each geometrical distance 
and the specification of the on-site terms jjl and ( [60, 61]. On the other hand, in 
a VBC state \xij\ 2 + |Ajj| 2 may be different from bond to bond and, therefore, the 
specification of VBCs amounts to giving the pattern of amplitudes of Xij an d Ay 
at each geometrical distance, in addition to specifying the SU(2) fluxes through the 
plaquettes. These parameters are the ansatze of a given state and serve as the variational 
parameters in the physical wave function that are optimized within the variational 
Monte Carlo scheme to find the energetically best state. It is worth mentioning that we 
use a sophisticated implementation of the stochastic reconfiguration (SR) optimization 
method which allows us to obtain an extremely accurate determination of variational 
parameters [62, 63]. Indeed, small energy differences are effectively computed by using a 
correlated sampling, which makes it possible to strongly reduce statistical fluctuations. 
This feature is especially important for the spin-1/2 QHAF since the energies of all the 
competing phases are rather close. 
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Figure 2: A hierarchical flowchart sorting out the myriad of different 6-, 12- and 36-site 
unit cell VBCs in order of increasing (from top to bottom) number of broken point group 
(PG) symmetry elements. The square boxes contain the VBC names followed by their 
respective symmetry PG. The 'parent' (maximally symmetric) VBCs are marked in red 
and those which have been found as competing ground states in studies using quantum 
dimer models are marked in pink [49, 50, 55]. The corresponding VBC patterns, and 
their discussion, are given in the text. As much as possible, we use labeling for the 
VBCs which is similar to that used in [55]. 



1.2. Parent spin liquid states 

The ansatz for the energetically best variational state, the U(l) Dirac SL, is given in 
figure la. Due to the U(l) flux if being and 7r [exp (i(p) = ripiaquctte Xij] through 
triangles and hexagons, respectively, it is denoted as [0, it] SL. In its mean-field band 
structure the Fermi surface collapses to two points at which the spectrum becomes 
relativistic with Dirac conical excitations [32]. Another energetically competing state, 
the uniform RVB SL, has zero flux through all plaquettes and is therefore denoted 
as [0, 0] SL. Its mean-field band structure consists of large circular spinon Fermi 
surfaces [34]. Both these states are fully symmetric, U(l) gapless SLs and can be 
extended to include second NN hoppings, leading to a gain in energy without changing 



Valence- bond crystals in the kagome spin-1/2 Heisenberg antiferromagnet 




9 



v. 



i /TV i 

i A A i 




(a) SVBC 



(b) VBC 



V 
(c) SVBC" 



Figure 3: The most symmetric 12-site unit cell VBCs: the center of symmetry is marked 
'C (the center of the shaded hexagon), around which bonds connected by the given PG 
symmetry operations are marked with the same color and style of the line. We will 
henceforth refer to these bonds as being in the same class, (a) The Star-VBC has the 
maximal PG symmetry, Cq v ; hence it acts as a 'parent' VBC. Its bonds breakup into 
three distinct classes, (b) The VBCo lacks crystallographic axes reflection symmetries 
in contrast to the SVBC; thus its symmetry group is reduced to C 6 . It has four classes 
of bonds, (c) The Star-VBC a has reduced (27r/3) rotation symmetry but preserves 
reflection symmetry; thus its symmetry group is C3,,. It has six classes of bonds. 



their nature [35]. It is worth noting that the effect of projection on these mean-field 
states can be drastic. 



2. Symmetry classification and enumeration of valence-bond crystals 

(VBCs) 

The VBC states on the kagome lattice break its elementary (three-site) unit cell 
translation symmetry with different unit cell sizes, which describe their modulation. 
In previous studies [51, 52, 53, 54, 43, 44, 45, 47, 49, 50, 55], using different methods, 
VBCs with 6-, 12-, 18- and 36-site unit cells were identified as possible ground states of 
the spin-1/2 QHAF. In this work, we will restrict our analysis to VBCs with 6-, 12- and 
36-site unit cells. For each unit cell size with a given center of symmetry, we enumerate 
VBCs starting from the maximally symmetric (Cq v ) 'parent' VBC and systematically 
break point group symmetry elements, right down to the VBC with no symmetry at all. 
This results in an enumeration of 19 VBCs in total, 9 VBCs each for 12- and 36-site 
unit cells and 1 VBC for the 6-site unit cell (see figure 2). Only 6 out of the 19 VBC 
have been studied previously. In this paper, we will study the possibility of any of these 
VBCs to occur as the ground state. 
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(c) VBC 3 




(d) VBC 



Figure 4: Other 12-site unit cell VBCs: (a) the VBCq has only reduced rotation 
symmetry (27r/3); thus in contrast to SVBC" its symmetry group is reduced to C3. 
It has eight classes of bonds, (b) The diamond- VBC has two perpendicular axes of 
reflection symmetry, thus giving rise to C^ symmetry, with seven classes of bonds, (c) 
The VBC3 has only 7T- rotation symmetry; thus its symmetry group is C2. It has 12 
classes of bonds, (d) The VBCi possesses only a single axis of reflection symmetry 
which bisects the sides of the shaded hexagon; consequently, its symmetry group is 
Cit,. It has 14 classes of bonds, (e) The VBC^ has the same symmetry as VBCi, but 
its reflection symmetry axis passes through a vertex of the shaded hexagon; we shall 
denote the symmetry group as C' lv to distinguish it from that of VBCi. It has 12 classes 
of bonds, (f) The VBC2 has no symmetry whatsoever; hence its symmetry group is just 
the identity, denoted here as E. Consequently, it has 24 distinct classes of bonds. 



2.1. 12-site unit cell VBCs 

The kagome lattice can be viewed as a triangular lattice of 12-site blocks shaped in the 
form of 'stars'. Within this picture, it was argued in [53, 54] that the ground state of the 
spin-1/2 QHAF has possible long-range singlet order that settles in this triangular star 
arrangement. This lends support to the picture that the ground state can be a VBC 
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(a) HVBC 



(b) HVBC 



(c) HVBC" 



Figure 5: The 36-site unit cell VBCs: the center of symmetry is marked 'C (the center 
of the shaded hexagon). The perfect hexagons, marked at their centers by 'P', form a 
honeycomb lattice at the center of which lie the shaded hexagons, f (a) The hexagonal- 
VBC has the maximal PG symmetry, Cq v ; hence it acts as a 'parent' VBC. Its bonds 
breakup into seven distinct classes, (b) The hexagonal-VBC , in contrast to the HVBC, 
lacks reflection symmetries about crystallographic axes; thus its symmetry group is 
reduced to Cg. It has 12 classes of bonds, (c) The HVBC Q has reduced (2ir/3) rotation 
symmetry but preserves reflection symmetry; thus its symmetry group is C3,,. It has 14 
classes of bonds. 



with a 12-site unit cell capturing some modulation. In total, nine symmetry distinct 
VBCs with a 12-site unit cell can occur; see figures 3 and 4 for their NN patterns. In 
particular, the SVBC state (figure 3a) was argued in [52] to occur as an instability of the 
U(l) Dirac SL and to be consequently stabilized as the ground state of the NN spin-1/2 
QHAF. Numerical studies using projected wave functions have shown this proposal to be 
incorrect and have also established the stability of the uniform RVB SL w.r.t. dimerizing 
into the SVBC state [32, 34, 35]. Furthermore, a recent QDM study [55] found the VBC 3 
state (figure 4c) to be a competing ground state and a DMRG study [41] concluded that 
the DVBC state (figure 4b) is close by in a generalized parameter space. In section 3, we 
study the possibility of a ground state realization of VBC 3 and DVBC states numerically, 
within the framework of projected wave functions. In fact, we perform this study for all 
12-site unit cell VBCs. 



I The points P are not centers of inversion (-7r-rotation) symmetry, as has been mismarked in figure 1(a) 
of [55], which corresponds to the HVBCq state in the present work. 
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(a) HVBC£ 



(b) 36-DVBC 




(c) 36-VBC 3 



Figure 6: (a) The HVBCq has only a reduced rotation symmetry (27r/3); thus in contrast 
to HVBC" its symmetry group is reduced to C3. It has 24 classes of bonds, (b) The 
36-diamond-VBC has two perpendicular axes of reflection symmetry, thus giving rise to 
C21, symmetry, with 19 classes of bonds, (c) The 36-VBC3 has only ir rotation symmetry; 
thus its symmetry group is C 2 . It has 36 classes of bonds. 



2.2. 36-site unit cell VBCs 

The building blocks of the kagome lattice can take nontrivial forms such as a 2\/3 x 2\/3 
expansion of the elementary 3-site unit cell, thus giving rise to a tilted 36-site unit cell. 
It was shown in [43, 44] that such a construction maximizes the density of hexagons on 
which dimer resonances occur, thereby lowering the energy. This 36-site unit cell VBC 
was studied numerically using series expansion [45, 46] and MERA [47], which found 
it to be a good approximation to the ground state of the NN spin-1/2 QHAF. Similar 
conclusions were also obtained from a QDM study [49, 50, 55]. Motivated by these 
findings we classify all 36-site unit cell VBC patterns on the kagome lattice, which leads 
to the identification of nine symmetry distinct VBCs; see figures 5, 6 and 7 for their NN 
patterns. 

In our previous work [35] we studied the HVBC state (see figure 5a) by using 
projected wave functions and found it to be higher in energy compared to the gapless 
SLs. However, the symmetry of the VBC identified in QDM studies [49, 50, 55] is that 
of the HVBCo state (see figure 5b), which has a lower symmetry compared to the HVBC 
state. In section 3, we study the possibility of a ground state realization of the HVBC 
state for the NN and NNN spin-1/2 QHAF. 
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(a) 36-VBCi 



(b) 36-VBC'! 



Figure 7: (a) The 36-VBCi possesses only a single axis of reflection symmetry, which 
bisects the sides of the shaded hexagon; consequently, its symmetry group is C\ v . It 
has 38 classes of bonds, (b) The 36-VBC^ has the same symmetry as 36-VBCi, but its 
reflection symmetry axis passes through a vertex of the shaded hexagon; we shall denote 
its symmetry group as C^, to distinguish it from that of 36-VBCi. It has 36 classes of 
bonds. Note that the 36-VBC2 (which has no symmetry at all) has not been drawn. Its 
symmetry group is just the identity E. Consequently, it has 72 distinct classes of bonds. 



2. 3. General remarks on the VBC classification 

It is worth mentioning that this VBC classification (for a given unit cell) is based 
on very general considerations of symmetry only and hence is not dependent on the 
formalism in which one studies these phases. In principle, it is possible to translate its 
construction from one language (e.g. QDM) into another (e.g. Schwinger fermions or 
bosons) for a VBC with a given symmetry. Moreover, within a given framework there 
can be different ways of constructing wave functions for a given VBC, consistent with its 
symmetry group. Firstly, one can add amplitudes beyond NN, consistent with the VBC 
symmetry group. Since we will study these phases within a slave particle approach, one 
can construct at the naive level simple mean-field wave functions or go much beyond 
mean-field and include the effects of full projection. At a next level, it is possible to 
improve the wave function by applying the Hamiltonian operator on it a given number of 
times and considering an optimized linear superposition of these wave functions with the 
original projected wave function. It is also worth noting that this hierarchical sorting 
of VBCs in each fixed symmetry sector also greatly eases the numerical search for a 
possible VBC stabilization as the ground state of the spin-1/2 QHAF. 
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3. Numerical Results 

We study the energetics of SL and VBC phases for the spin-1/2 QHAF using Gutzwiller 
projected fermionic wave functions with the variational quantum Monte Carlo technique. 
Our variational calculations are performed on clusters with 432 (i.e. 3 x 12 x 12) or 
576 (i.e. 36 x 4 x 4) sites and mixed periodic-antiperiodic boundary conditions which 
ensure non-degenerate mean-field wave functions at half filling. The large size of the 
cluster ensures that the spatial modulations induced in the observables by breaking 
of rotational symmetry (due to mixed boundary conditions) remain smaller than the 
uncertainty in the Monte Carlo simulations. 

Among the class of NN fully symmetric and gapless SLs, the U(l) Dirac SL ([0,7r] 
SL) has the lowest energy for the NN spin-1/2 QHAF. On a 432-site cluster its energy per 
site is E/Ji = -0.428 63(2) and the uniform RVB SL ([0,0] SL) has a slightly higher 
energy per site, E/J\ = —0.412 16(1) [32]. For the 576-site cluster these values are 
E/J-y = -0.428 66(1) for the U(l) Dirac SL and E / J x = -0.411 97(1) for the uniform 
RVB SL [35]. Upon inclusion of NNN hopping amplitudes, one gets the extended U(l) 
Dirac SL or the extended uniform RVB SL, which are labeled by one additional flux 
through a plaquette of the type '234' in figure la, the flux through the other triangular 
plaquette formed by NNN bonds only is then fixed. Hence, the extended Dirac SL can 
be either the [0, 7r; tt, 0] or the [0, 7r; 0, %] SL and analogously the extended uniform RVB 
SL can be either the [0, 0; tt, tt] or the [0,0; 0,0] SL [35]. For the NN spin-1/2 QHAF 
these extended SLs have a slightly lower energy, but they perform much better for the 
J\ — J 2 spin-1/2 QHAF, see figure 9. 

3.1. Results on the stability of gapless SLs towards VBC perturbations 

We carried out an extensive numerical study of the local and global stability of the 
NN U(l) Dirac and uniform RVB SL toward dimerizing into all 6-, 12- and 36-site 
unit cell VBCs. In cases where we did find dimerization with NN bond amplitudes, 
we added second NN bond amplitudes to the SL and VBC ansatz (consistent with 
symmetries), since this led to a significant gain in energy. Our main focus was on the 
CVBC (figure lb), DVBC (figure 4b), VBC 3 (figure 4c) and HVBC (figure 5b) states, 
since these have been identified as ground states of the spin-1/2 QHAF in other studies. 
We perform our analysis by first fixing a background flux corresponding to the SL liquid 
whose stability we wish to study. Then, we introduce an amplitude modulation of x%j 
consistent with the symmetries of the VBC, i.e. bonds belonging to the same class 
(color/line marking in figures lb, 3, 4, 5, 6 and 7) have the same amplitude (xx)i which 
is set to different values for different classes. Starting from an arbitrary unbiased point 
(Xa's) in the variational space we perform an optimization of the wave function to obtain 
the lowest energy state [62, 63]. 

3.1.1. The case of the U(l) Dirac SL For the NN spin-1/2 QHAF, the variation of 
parameters and energy in the Monte Carlo optimization for the four competing VBCs 
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Figure 8: A typical variational Monte Carlo optimization run for the CVBC (a), DVBC 
(b), VBC 3 (c) and HVBC (d) wave functions, for the NN S = 1/2 QHAF. The 
variational parameters xx an d energy (insets) are shown as a function of Monte Carlo 
iterations. The NN U(l) Dirac SL corresponds to \xx\ = 1- On starting from different 
sets of initialized parameter values we return back (within error bars) to the U(l) SL. 
The optimized parameter values are obtained by averaging over a much larger number 
of converged Monte Carlo iterations than that shown above. 



(regarded as a dimerization of the U(l) Dirac SL) mentioned above is given in figure 8. 
As can be clearly seen, the energy converges neatly to the reference value of the U(l) 
Dirac SL, and all the parameters converge to xx — 1 (within error bars) after averaging 
over a sufficient number of converged Monte Carlo steps; thus the translation symmetry 
associated with the SL is restored. In fact, we performed these calculations for all 
6-, 12- and 36-site VBCs and found that in each case the U(l) Dirac SL is stable 
towards opening a gap and destabilizing into any of these VBCs. This remarkable 
stability (for all VBCs) is also preserved upon addition of a NNN ( J2) super-exchange 



Valence- bond crystals in the kagome spin-1/2 Heisenberg antiferromagnet 



13 



-U.^ 


11 




1 1 1 1, 

L 


»*'- 

""■"Oow-.. : 




^ -0.41 

"03 

<3 -0.42 
>> 




[0,7i;ji,0] 

O O [0,7i;0,Ji] 

A- - -A [0,0;7C,Jl] 

A A [0,0;0,0] 

□ □ HVBC () 

[0,7C] 

[0,0] 


ISJ 


~ ~'y^' 


<B -0.43 
-0.44 





J, =-0.045 

2,c 



-0.3 -0.2 -0.1 A 0.1 

L (FM) L (AF) 



Figure 9: Energy versus J 2 for SLs and the HVBC state (see figure 5b). The HVBC 
state becomes the lowest in energy for J 2 ~ —0.045. Error bars are smaller than the 
symbol sizes. 

coupling in the Hamiltonian of both AF and FM type. We verified these results by 
doing many optimization runs starting from different initial values of the parameters in 
the respective variational spaces. Thus, we can safely conclude that the U(l) Dirac SL 
has the lowest variational energy among all proposed competing VBC states, at least 
within the Schwinger fermion representation of the spin model for J 2 greater than a 
certain critical value J 2c , which is given and discussed in the ensuing text. 



3.1.2. The case of the uniform RVB spin liquid We now shift our focus to the uniform 
RVB SL and address the question of its stability. For the NN and NNN (AF and FM) 
spin-1/2 QHAF, we find that all 6- and 12-site unit cell VBCs have a higher energy 
compared to the uniform RVB SL. However, interestingly enough, for the NN spin-1/2 
QHAF, this NN uniform RVB SL opens up a gap and destabilizes into a 36-site unit cell 
VBC, namely the HVBCo state (see figure 5b). The gain in energy due to dimerization 
becomes more pronounced on addition of second NN hopping amplitudes to the wave 
function which are consistent with Cq symmetry. On adding a NNN super-exchange 
coupling of FM type to the Hamiltonian and following this second NN HVBCo state 
(now, a dimerization of the extended uniform RVB SL), one finds that it becomes the 
lowest in energy for J 2 ~ —0.045 (see point A in figure 9), consistent with the findings 
in [59]. It is worth noting that the symmetry of this VBC is precisely that of the VBC 
identified in the QDM study [49, 50, 55] and has a lower symmetry compared to the 
HVBC state that was previously studied by us with similar conclusions [35]. The flux 
pattern of this VBC consists of flux through all elementary triangles, hexagons and a 
71 flux through the '234 plaquettes (see figure la) inside the perfect hexagons only. The 
lower symmetry of the HVBCo compared to the HVBC implies a larger variational space 
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of hopping amplitudes and consequently a lower energy that is seen from the fact that 
the level crossing or the onset of VBC order is shifted from J 2 ~ —0.09 [35] for HVBC 
to J2, c ~ —0.045 for HVBCo state. Thus our results still point to a gapless ground state 
for J2 ~ —0.045, which is along the lines of our previous work [35, 36]. 

4. Conclusions and discussions 

In this paper, we enumerated all 6-, 12- and 36-site unit cell VBCs based on symmetry 
considerations alone and subsequently investigated the possibility of stabilizing any of 
these VBCs in the NN and NNN spin-1/2 QHAF on a kagome lattice. We found that 
the U(l) Dirac SL is remarkably robust toward dimerizing into any of these VBCs, for 
both the NN and NNN spin-1/2 QHAF. However, the uniform RVB SL dimerizes into a 
36-site unit cell VBC, which becomes the lowest in energy on addition of a very weak FM 
coupling, J2 )C ~ —0.045. Our systematic and thorough numerical investigation brings us 
to the conclusion that, at least within the Schwinger fermion approach to the spin model, 
the U(l) Dirac SL has the best variational energy for J 2 ^ —0.045. The conflict between 
our results, which point to a gapless ground state in this region, and those obtained by 
exact diagonalizations and DMRG calculations, which instead suggested the presence of 
a fully gapped spectrum, remains open and deserves further investigation. One possible 
direction would be to include vison dynamics in the projected wave functions [64], 
which may be necessary to capture topological order faithfully. Another step would be 
to improve our variational wave functions based on the application of a few Lanczos 
steps [65] and then perform an approximate fixed-node projection technique. The 
possibility that an unconventional VBC breaking time-reversal symmetry is stabilized 
as the ground state cannot be ruled out [55]. Finally, we mention that VBC order might 
also set in via confinement transitions of the Z2 SLs [66] , this remains to be investigated 
numerically. 
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